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We present the first attempt to use a combination of CMB, LIGO, and PPTA data to constrain 
both the tilt and the rnnning of primordial tensor power spectrnm through constraints on the 
gravitational wave energy density generated in the early universe. Combining measurements at 
different cosmological scales highlights how complementary data can be used to test the predictions 
of early universe models including the inflationary consistency relation. Current data prefers a 
slightly positive tilt (n^ = O.Obi'lQ gg) and a negative running {nt,iun < —0.22) for the tensor power 
spectrum spectrum. Interestingly, the addition of direct gravitational wave detector data alone puts 
strong bounds on the tensor-to-scalar ratio r < 0.2 since the large positive tensor tilt preferred 
by the Planck temperature power spectrum is no longer allowed. Adding the recently released 
BICEP2/KECK and Planck 353 GHz polarization cross correlation data gives an even stronger 
bound r < 0.1. We comment on possible effects of a large positive tilt on the background expansion 
and show that depending on the assumptions regarding the UV cutoff (fcuv/fc* = 10^^) of the 
primordial spectrum of gravitational waves, the strongest bounds on nt = 0.07to'go derived 
from this effect. 

PACS numbers: 


I. INTRODUCTION 

The detection of B-mode polarization on large angu¬ 
lar scales by the Background Imaging of Cosmic Extra- 
galactic Polarization II [U BICEP2] collaboration, and 
its possible primordial origin as the result of relic gravita¬ 
tional waves has invigorated the cosmological community. 
There has been much debate about the interpretation of 
the BICEP2 data, and this is an important issue that will 
be resolved with independent measurements and multi¬ 
frequency observations. The recent BICEP2/KECK and 
Planck 353GHz polarized map joint analysis has shown 
that at least a part of the BICEP2 and KECK signal is 
due to polarized dust [2]. However, more importantly the 
observations have yielded renewed interest in constrain¬ 
ing models of the early universe through their predictions 
of the spectrum of primordial gravitational waves and the 
testability of the inflationary scenario. 

Fortunately, at least for a simple class inflationary 
models in which a single scalar field drove a period of 
inflation more or less directly preceding the radiation 
dominated era of hot big bang expansion of our uni¬ 
verse, positive predictions can be used to test the in¬ 
flationary paradigm. The primordial power spectrum of 
tensor modes is typically parameterized as a power law 
0 : 

pm = A{k.) , ( 1 ) 

where At is the tensor amplitude, k is the wavenumber, 
fc* is some reference wavenumber, and nt is the tensor 
spectral index, or tilt. While there are, in general, de¬ 
viations from this simple expression, it serves to param¬ 
eterize a host of scenarios. For the simplest single field 
models of inflation there is a relation between the tensor 


spectral tilt and the tensor-to-scalar ratio r, sometimes 
referred to as the inflationary consistency condition [3]. 

At leading order in slow-roll, this relation is given by 

[51 [S] 

r=^= -8nt , (2) 

where is the amplitude of the scalar power spectrum. 
In more complicated models with, for example, a non¬ 
trivial sound speed or multiple dynamical fields, this re¬ 
lation is modified to be r < —Sn* (see e.g. [ZHS] for 
recent discussions of some of these issues). In all mod¬ 
els of inflation where the tensor perturbations originate 
as vacuum fluctuations, the tensor amplitude At is set 
only by the expansion rate which is nearly constant dur¬ 
ing inflation and always decreases when the null energy 
condition holds: inflation always predicts a small neg¬ 
ative tensor tilt. The specific value of the tilt nt may 
depend on the details of the model, but a negative ten¬ 
sor tilt is a generic prediction of inflation when the null 
energy condition holds. There are however, other early 
universe models which make different predictions for the 
tensor spectrum [IDHIl]. Hence we take a phenomeno¬ 
logical point of view and do not restrict our analysis to 
negative values of nt- If the data provide evidence for a 
positive tensor tilt, we should further explore alternatives 
to inflation. A detection of a positive tensor tilt would 
suggest at the very least that the primordial gravitational 
waves were not produced as vacuum fluctuations during 
inflation^. 


^ That being said, counter examples do exist even in the context of 
inflation 11311151 . but fortunately these models tend to make addi- 
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The aim of this paper is to perform the first joint anal¬ 
ysis of cosmic microwave background (CMB), Laser In¬ 
terferometer Gravitational wave Observatory [151 LIGO], 
and Parkes Pulsar Timing Array [T71 PPTA] data to put 
a constraint on both the tilt nt and and the running 
n-t,run of the primordial tensor power spectrum. Previous 
studies have shown that such combinations of data would 
improve the available lever arm of the data [TH] ■ BICEP2 
data alone yields only a weak constraint on the tensor tilt 
nt as the polarization signal of gravitational waves drops 
rapidly as a function of multipole; hence the CMB alone 
provides a limited lever arm to constrain the tensor power 
spectrum. We will show that the constraint on nt is dom¬ 
inated by the Planck temperature power spectrum data, 
which prefers lower power on large angular scales^, result¬ 
ing in a preference for positive tensor tilt. The additional 
data extend our lever arm significantly by increasing the 
scales over which we can measure the gravitational wave 
energy density by 18-20 orders of magnitude. Although 
the constraints at the small scales probed by LIGO and 
PPTA are only upper limits on the gravitational wave 
energy density, earlier estimates have shown [251 HI] that 
these experiments are particularly sensitive to the tilt of 
the primordial tensor mode power spectrum, assuming 
a detection of a tensor-to-scalar ratio r > 0.01. The re¬ 
sults presented in this paper are relevant for constraining 
the details of the gravitational wave power spectrum, in¬ 
dependent of the current interpretation of the BICEP2 
data. We will perform each analysis with and without 
BICEP2 data, and show that one can get meaningful re¬ 
sults even without constraints from BICEP2 on r. In 
addition, we will include the recent released cross anal¬ 
ysis data of BICEP2/KECK and the 353 GHz Planck 
polarization data [2|. 


The paper is organized as follows. In Sec.jTTjwe briefly 
review the relevant physics. In Sec. [ml we present our 
results on the bounds on nt and from a variety of 

data combinations, including in Sec. El a discussion of 
the observational constraints derived from the impact of 
the primordial gravitational wave energy density on the 
background expansion of the universe not included in the 
initial analysis. We summarize some future prospects and 
conclude in Sec. |IV| A summary of the derived parameter 
constraints using various combination of data sets can be 
found in Tab. [I] All bounds quoted in the paper are 95% 
bounds. A discussion of the definition of gravitational 
wave energy density is included in Appendix [A| 


tional predictions for other observables such as non-Gaussianity 
which would allow additional checks of these scenarios. 

^ To put any reasonable constraints one has to take a pivot scale 
that is within the observable window. In this paper we chose the 
slightly unconventional /c* = 0.01 Mpc~^. See e.g. Ref. |19| for 
a discussion of pivot scale in the context of constraining gravita¬ 
tional waves with B-mode polarization. 


II. GRAVITATIONAL WAVE SPECTRUM AND 
ENERGY DENSITY 

While measurements of the B-mode polarization of the 
CMB constrain the tensor mode power spectrum directly, 
observations from gravitational wave observatories and 
pulsar timing arrays instead place limits on the energy 
density of gravitational waves at specific wavenumbers 

mm 
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Ptjk) 
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Tir]o,k) , 


( 3 ) 


where rjo is the conformal time in the current epoch, 
is the tensor transfer function, and an over¬ 
dot refers to a derivative with respect to cosmic time 
t. The transfer function projects the primordial signal 
from early to late times and is obtained by solving the 
Klein-Gordon equation for the tensor fluctuations in an 
expanding universe^. See Appendix]^ for a derivation of 
this formula. 

In the case of single field slow-roll inflation, the power 
spectrum of the tensor modes, Pt{k), is well described 
by Eq. 0 with ut ~ 0. The leading deviation from 
the simple power law Eq. 0 expected in some models is 
parameterized by allowing the tilt to vary with scale: 


/ h \ "t(fc.)-|-5nt,run ln(fc/fe,) 

Pt{k) = At{K)lyj-] 


( 4 ) 


where nt^mn = dnt/dlnk is the running of the tensor 
spectral index and again At is the amplitude of the ten¬ 
sor power spectrum. Single field slow-roll inflation makes 
a model-dependent prediction for the tensor-to-scalar ra¬ 
tio, i.e. r = At I As = 16e [3]. The (single field slow- 
roll) inflationary consistency condition also tells us that 
nt = —r/8, and similar relations can be derived for the 
running and higher order terms which are generally slow- 
roll suppressed [25] . As a result, the scale dependence 
of tensor fluctuations is described by a limited number 
of parameters, and should be uniquely determined by a 
limited number of measurements. Despite being slow-roll 
suppressed, higher order terms can become quite signifi¬ 
cant in determining the spectrum far from the pivot scale 
fc* |26] . We are not yet in the data-driven regime, how¬ 
ever, and with only upper limits from PPTA and LIGO 
on the gravitational wave energy density, the data is not 
yet capable of falsifying the consistency conditions with 
high significance. Hence, we include the running but 
leave a treatment of higher order corrections to future 
work. 

In this paper we assume that the power law form of 
the tensor power spectrum holds all the way up to some 


® In general a non-standard thermal history will affect the form of 
the transfer function. See for example 12111241 . 
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fixed ultraviolet cutoff, which we will discuss in more de¬ 
tail below. Physically, this means that we are assuming 
an instantaneous transition from the phase of the early 
universe responsible for producing the primordial spec¬ 
trum to the phase of radiation domination. This instant 
reheating assumption is likely to be violated in a realistic 
model of the early universe where one would expect modi¬ 
fications to the power spectrum at wavenumbers near the 
cutoff (see [ 21 ] for an excellent review). Other effects, 
such as the late production of entropy can also modify 
the spectrum and weaken the constraints derived below. 
Hence our constraints apply to a phenomenological model 
with the simplest assumptions: a simple power law pri¬ 
mordial spectrum (including the possibility of running), 
instantaneous reheating, and a standard thermal history. 
Any specific model which violates one or more of these 
features is likely to be subject to different observational 
constraints. 

We use a pivot scale of fc* = 0.01 Mpc“^, particularly 
because the B1CEP2 experiment is insensitive to smaller 
scales [12]. We use an analytic solution for the gravita¬ 
tional wave transfer function valid at late times, which 
has an accuracy of about 1 % in a flat universe and is 
relevant for very small scales k ^ fcgg and late times 
‘n Veq m- We modified the publicly available Markov 
Chain Monte Carlo sampler cosmomc |31| and the Boltz¬ 
mann code CAMS [32], by adding a module^ that com¬ 
putes the gravitational wave energy density as a function 
of wavenumber k, which can be converted to frequency 
using the relationship 


^ 27ra(?7o) ’ 


The observations of pulsar timing arrays correspond 
to frequencies of / €(10“® Hz, 10“^ Hz), while 

gravitational wave observatories probe / g( 0.1 mHz, 
5000 Hz) [ 33 ]. The constraint from LIGO^ is ffcw < 
5.6 X 10“® ( 68 /Ho)^ [ 31 ] at 2cr, while the upper limit from 
PPTA is flaw < 1-5 X 10"® (68/Ho)2 [ 33 ] at 2cr where 
Hq = 100 hkms“^Mpc“^. 

For each value of the gravitational wave energy den¬ 
sity, a likelihood is computed based on the upper limits 
from PPTA and LIGO (y^ = 4.0, corresponding to 2a 
at 95% confidence level) from a truncated, logarithmic 
distribution. We use an analytic solution to obtain the 
conformal time 77 as a function of scale factor [36] . 


^ Code will be made publicly available through http; //www. astro. 
princeton.edu/'meerburg/coding/ 

We have approximated the constraints from LIGO by taking a 
constant upper limit over the frequency range considered. For 
the purpose of this analysis, this approximation is sufficient. 



FIG. 1: The marginalized 2D posterior constraints on nt and 
r. The large grey contonr shows the constraints of using BI- 
CEP2 in conjunction with CMB data from Planck [28], the 
HST prior on the Hubble constant [29], and the constraints 
on the aconstic scale from BAO |30| measurements (the lat¬ 
ter three combined in the name ‘All’). Adding in data from 
the PPTA constraints on the gravitational wave energy den¬ 
sity improves the constraints three-fold on nt, while the con¬ 
straint on r is largely unchanged. Using data from LIGO 
instead tightens this even further, highlighting the constrain¬ 
ing power of LIGO on the primordial tensor power spectrum. 
All contours agree on the value of r announced by the BI- 
GEP2 team, confirming that the LIGO and PPTA data are 
not particularly sensitive to the value of r directly, but they 
are sensitive to the tilt nt, as shown in Eq. (|^. 


III. RESULTS 

A. Combining CMB measurements with 
low-redshift constraints 

In all runs we include temperature data from the 
Planck satellite |28j with the polarization prior on the 
optical depth from WMAP 127 ], WMAP low £ polariza¬ 
tion, a prior on the Hubble parameter from the HST key 
project [22] and BAO measurements from the Sloan Dig¬ 
ital Sky Survey [SO] SDSS]. We then include a variety 
of other probes, including the recent measurements of 
large-scale polarization from the BICEP2 team [T], and 
constraints on flew from the PPTA or from LIGO data. 
In order to be insensitive to lensing effects in polarization 
data, we restrict the BIGEP2 data to the four lowest mul¬ 
tipole bins. Last, we include an analysis combining the 
previous data sets with the BIGEP2/KEGK and Planck 
cross correlation likelihood [ 2 ] . Here we fit for the lensing 
amplitude and the dust contribution. 

Our results are presented in Fig. [^ and Table |l] As ex¬ 
pected, we get a weak constraint on the 95% limit value 
of nt of 1.64 < nt < 2.63 when using only BIGEP2 data. 
The constraint on nt is improved signihcantly by adding 
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Dataset 

Parameter 


nt 

ro.oi 

Tit,run 


Best fit 

Mean, 95% limits 

Best fit 

95% limits 

Upper, 95% limits 

All 

2.66 

2.08^;j® 

0.24 

< 1.0 

- 

All+LIGO 

0.51 

0.i3«:®4 

0.06 

< 0.2 

- 

All+PPTA-bLIGO (running) 

2.99 


1.06 

< 1.7 

< -0.32 

All {N,s,kvv/k, = 10^°) 

0.52 

0.12«:« 

0.01 

< 0.21 

- 

All {N,a,kvy/k, = 10^4) 

0.44 

0.05l°;f 

0.01 

< 0.17 

- 

AlH-BKxP+PPTA 

0.65 

0.32+gJ^ 

0.04 

< 0.1 

- 

AlH-BKxP-tLIGO 

0.26 

0.06t°:“ 

0.02 

< 0.1 

- 

All-fBKxP-tPPTA+LIGO (running) 

3.96 


0.02 

< 0.11 

< -0.22 

All-bBKxP (A'eff,A:uv/fc. = 10^°) 

0.49 

n 1 c+O.SS 

U.iO_o.80 

0.05 

< 0.09 

- 

All+BKxP (Aeff,/fcuv/fe* = 10^^) 

0.39 

0 n7+0-^2 

U.UI o.80 

0.04 

< 0.09 

- 

All-tBICEP2 

1.42 

-l-*OU_i 03 

0.18 

o.i8«:J9 

- 

All-tBIGEP2-tPPTA 

0.66 

n rt::+0.29 
U.OO_o.4g 

0.22 

f) 1 q+0.09 

u.iy_o.o8 

- 

All-tBIGEP2-tLIGO 

0.56 

0 ‘^p;+0-20 

U.OO_o.3g 

0.18 


- 

All+BICEP2+PPTA-bLIGO (running) 

0.56 

1 1 r+i.ae 

21 

0.19 

0.19l°;l 

< -0.05 

All-tBIGEP2 (iVeff,fcuv/fc, = 10^°) 

0.47 

n 4+0 09 

^•+-0.35 

0.187 

f) 1 7+0.07 
' -0.06 

- 

All-tBIGEP2 {N,b, feuv/fc* = lO^'^) 

0.41 

n OO+ 0.12 

u.oo_0.25 

0.16 

n 17+0.07 
< -0.06 

- 


TABLE I: Summary of parameter constraints using various combinations of data sets. “All” refers to Planck (2013), WMAP 
low i polarization, HST and BAO data. In all studied cases a positive nt is preferred over a negative n^, but the significance 
is small when you remove any BICEP2/KECK data. Recall that we use a pivot scale fc* = 0.01 Mpc~^. 



'' 0.01 

FIG. 2: The marginalized contours in the r — nt plane com¬ 
paring BICEP2 to the BKP data. If all the signal in BICEP2 
is primordial, we obtain strong bounds on nt, with a 2a pref¬ 
erence for nt > 0. The upper limit is driven by the LIGO 
constraint. From the BIGEP2/KEGK cross Planck analysis r 
is lowered and therefore in combination with the constraints 
from Planck TT power spectrum leads to a shift in the peak 
value of nt- With a small value of r, constraints on nt are 
weakened. However, including LIGO puts a strong upper 
bound on nt as before. 


LIGO, lowering the 95% upper limit by roughly a fac¬ 
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FIG. 3: The marginalized contours in the r — nt plane without 
using BIGEP2 data. When ro.oi = 0 there is no constraint on 
nt in principle, however Planck TT can be best fit with ro.oi = 
0.24, hence there exists some constraint on nt- Planck prefers 
a very blue tensor spectrum, in order to create an artificial 
running of the temperature data on large scales. From this 
analysis it is also apparent that Planck drives the value of 
Tit, not BIGEP2 data. Interestingly, by adding LIGO, such 
large values of nt are ruled out, which has the effect that r is 
significantly constrained ro.oi < 0.2 at 95%. 


tor of 5 from 2.63 to 0.59. The mean value we get for 
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the tensor tilt at CL 95% is n* = 0.394lg;2g8. The up¬ 
per limit is a very hard cut, because of the exponential 
drop off. Our results reflect the knowledge that it is very 
difficult to measure n* from CMB data alone (see e.g. 
Refs. [HIST]) and small scale GW experiments are nec¬ 
essary to test the inflationary consistency condition. We 
are also aware that a detection of gravitational waves 
by LIGO will not be particularly useful for determining 
the primordial tensor spectrum as the signal on scales to 
which LIGO is sensitive is dominated by astrophysical 
sources of GW, such as binary mergers. Since current 
constraints are mainly upper limits on the gravitational 
wave energy density and any signal which is detected at 
LIGO is expected to be of non-primordial nature, what 
we present here are conservative constraints on Uj. 

When we replace the BICEP2 data with the cross data 
[2], there is a very strong upper bound on the tensor-to- 
scalar ratio from BIGEP2/KEGK data, as shown in Fig[^ 
We find rit = 0.06lg;8g and ro.oi < 0.1. 

Even after cross checking Planck with BIGEP2/KEGK 
data, there is still significant uncertainty in the total 
fraction due to dust. For completeness, we also consid¬ 
ered constraints using only Planck temperature data and 
WMAP polarization, as shown in Fig.[^ In principle, no 
bound can be put on nt when r = 0, but r = 0 is as likely 
as any other small value®. It is interesting to note that 
Planck drives the value of nt, not the BIGEP2 data. In 
fact, allowing nt to vary freely, leads to the best-fit value 
^ 0.01 = 0.24. We show the contributions from the best-fit 
spectra to the TT power spectrum at low multipoles in 
Fig. a For very positive tensor tilt, the TT data can 
be made more consistent with the theory by providing 
a possibility to lower the power on large scales. This is 
also the reason why Planck alone prefers a large value of 
r (to increase the overall contribution from the tensors 
to the TT spectrum) and with a large positive tensor 
tilt (to mimic suppression of the power on large scales). 
When LIGO data is added to the likelihood, such prefer¬ 
ence completely disappears, since nt can not be too large 
for most values of r. As a result, a very large value of 
r is no longer allowed and we actually get a constraint 
on r that is comparable to the constraint from r with¬ 
out adding the tilt and the constraints after adding the 
BIGEP2/KEGK and Planck cross data. As such, in the 
scenario that even a larger fraction in the BIGEP2 patch 
is due to dust, simply adding LIGO puts a very strong 
bound on r. That being said, obviously if all measured 
signal in BIGEP2/KEGK is due to dust, this puts a very 
strong bound on r. It is interesting to see that even with¬ 
out BIGEP2/KEGK data one could have put a constraint 
on r when allowing nt to vary; adding LIGO tells us that 


We do not encounter serious issues with convergence. An alter¬ 
native procedure is to constrain r at two different scales and to 
treat the tilt as a derived parameter. Note that additional scales 
would be required for constraints on the running and higher order 
scale dependence. 


r < 0.2, which is consistent with the BICEP2/KECK 
cross Planck data. 

Next, we consider the running of the tensor tilt. Be¬ 
cause the PPTA and LIGO measurements are at different 
scales we can use both PPTA and LIGO in our analy¬ 
sis. The results are presented in the right-hand panel of 
Figurej^ Without small scale constraints on the gravita¬ 
tional wave energy density data are incapable of putting 
any bound on the nt,run- As is shown in the right hand 
panel of Figure the running is exclusively confined 
to negative values and this constraint comes solely from 
PPTA and LIGO. Naturally, as nt increases, nt,run de¬ 
creases in order to compensate for the large values of the 
tilt. This causes the constraint on nt to worsen as one 
moves to lower values of the running. For values above 
nt - 3 (with BIGEP2) or nt - 6 (with BIGEP2/KECK 
cross Planck), the GMB measurements control the con¬ 
straint, and the joint data does not allow for larger values 
of the tilt. Since LIGO and PPTA do not provide lower 
limits, there is no useful constraints on the lower limit 
of the running. The constraints on the tensor-to-scalar 
ratio r are weakened, since running adds another free pa¬ 
rameter that allows the fit to improve on large scales in 
the TT power spectrum. 


B. Gravitational wave contribution to the massless 
degrees of freedom 


Our analysis has shown that the inclusion of small scale 
gravitational wave constraints results in a bound on nt ^ 
0.5 (at Icr). Since gravitational waves also contribute 
to the total radiation energy density of the universe, a 
consistent analysis must take into account their effect 
on the background expansion, as was considered in [2T] 
and [ 35 ] • The increased radiation energy density alters 
the peak structure of the GMB as explained in e.g. [ 53 ] 
and m and also affects the prediction for the primordial 
abundance of light elements, see for example [41] . Rossi 
et al. |42j recently presented constraints from Lyman- 
a forest measurements in conjunction with a variety of 
probes. 

The energy density of gravitational waves is given by 
(see Appendix [A] for details) 


PGw = 


327rGa2 Jr 


( 6 ) 


Here we used that the conformal time derivative of the 
transfer function during radiation domination is given 
by T'(r\ < 77eq,fc > fceq) = —kji{kr]) [33]. We cannot 
simply integrate over all wavenumbers, since this inte¬ 
gral diverges in the IR for nt < —4 and in the UV for 
nt > —2. Physically, one expects [44] that superhorizon 
modes should not contribute to the local energy density 
(note that even this itself is an ambiguous statement in 
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FIG. 4: The contribution from the tensors to the TT power spectrum. A very large positive tilt leads to a rescaling of the 
spectrum, except on very large scales, where the spectrum is suppressed. This suppression is what Planck prefers (though not 
significantly) over no tensor modes at all. 



FIG. 5: Gonstraints allowing for more freedom in the GW power spectrum. Left panel - the marginalized 2D constraints in the 
Ut — r for the cases including BICEP2 (bottom) and BICEP2/KEGK cross Planck (top) data. Right panel - constraints on 
the nt — nt,run for the same cases. Allowing for running in the tensor spectral index does not change the constraint on the tilt 
significantly, as the errors on the running are large and there is no signihcant correlation between the two parameters. BIGEP2 
data improves constraints on nt as expected. The addition of nt,run generally weakens the constraint on r. Note the different 
ticks in the figures on the right. 


the sense that the total energy density could only be mea¬ 
sured after averaging over several wavelengths |22] ) and 
this behavior is captured by the fact that small wavenum¬ 
bers make a negligible contribution to the integral. Fur¬ 
thermore, we do not expect that the physics of the early 
universe generated primordial gravitational waves at ar¬ 
bitrarily small scales. The largest conceivable wavenum¬ 


ber at which primordial gravitational waves could have 
been produced in the early universe is given by the Planck 
scale at the beginning of the hot big bang expansion of 
the universe, since modes with larger wavenumbers would 
have been super-Planckian at some point in the history of 
the universe. If we assume that the mechanism which is 
responsible for solving the flatness and horizon problems 
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FIG. 6: Left: The relative change in the scale of matter-radiation equality as a function of the tensor tilt rit, assuming r — 0.2. 
The simple analytical estimate suggests that nt ~ 0.4 would be modifying the scale of equality at the 5% level, which is close 
to current bounds from WMAP and Planck. Obviously, for smaller r, this constr aint wea kens accordingly. Right: The change 
in Ned as a function of nt for r = 0.2 using the two approximations given in Eqs. (111 and (12 I. The difference is very small up 
until the energy density in gravitational waves becomes of the order of the total energy density in the Universe. Our constraints 
are derived using Eq. ©• A detailed explanation of both curves and how these are derived can be found in the text. 


is also responsible for generating the primordial spectrum 
of gravitational waves, then we would expect that the 
spectrum should extend over a range of scales equiva¬ 
lent to the amount of expansion required to solve those 
problems. Below we give results for the case where the 
spectrum extends only to LIGO scales and also for the 
case where the UV scale is 60 e-folds smaller than the 
present horizon size. It is also worth mentioning that 
although the scaling of the energy for each mode k is 
equivalent to that of radiation, if there is some comoving 
IR bound, this scaling is (weakly) broken by the large 
scale modes which enter the horizon and begin to oscil¬ 
late. For our purposes, we are only interested in positive 
nt, for which the UV modes dominate the total energy 
density. It is natural then to introduce a new variable 
q = k/kuY, dk = kuydq and feip/Zcuv = £■ In the cases 
we consider, we can safely take the limit e —)■ 0, but our 
results will still depend on the choice of fcuvj which we 
leave unspecified for the moment. 

We can then perform the integral for positive nt 

( \ 

2^i^ + CI(l/(fcuv^))- (7) 

Hence for the gravitational wave energy density we find 
(up to corrections of order e) 


PGwiv) 


AsT J_^ 

327rG \ fc* / 2nt {ariY 


( 8 ) 


It is useful to relate this energy density to quantities 
which may be more familiar (for a similar analysis, 
see [45]). Deep in the radiation-dominated era, after 
big bang nucleosynthesis (BBN) but well before matter- 
radiation equality, total energy density of the universe is 


given by 


Ptot = P't + Pu + PGW 


— P~f 





(9) 


where = 3.046 in the standard model in the absence 
of gravitational waves |in|. Using the fact that during 
this period, l/(a? 7 )^ = H'^ = 87rGptot/3, and also assum¬ 
ing that the gravitational waves are not the dominant 
source of energy density, we can write the fraction of to¬ 
tal energy density in gravitational waves as 


Pgw _ AgV / fcxjv \ 

ptot 24nt V ^* / 

^ (iVeff- 3.046) p^ 

which can be solved for iVefi 


S (n)4/3 

AsC { 
24nt V 

T)”*' 

+ 3.046 

1 - 

AsP 

24nt 




( 10 ) 


( 11 ) 


This expression breaks down when the gravitational wave 
energy density becomes comparable to the total energy 
density of the universe. A fully consistent analysis for 
larger values of the gravitational wave energy density 
would require deriving a new transfer function which 
does not assume radiation domination, but instead in¬ 
corporates the backreaction of the gravitational waves 
on the expansion history. Luckily, the data is sufficient 
to constrain the gravitational wave energy density to a 
level significantly below the regime where this presents 



















a problem. An alternative approximation for A^eff which 
does not have a singularity can be found by neglecting 
the contribution of pcw to H in Eq. ([^ which gives the 
same result as the first order Taylor expansion of Eq. (11) 
about pGwIPtot = 0, 


Aeff « 3.046 + 3.046 



AsT /fcuvV 
24nt V fc* / 


( 12 ) 

Both Eq. ( |11[ ) and Eq. ( |12[ ) should only be trusted in the 
regime pcWptot < 1. We use Eq. 0 in our analysis 
since it gives more conservative constraints on rit , though 
the difference is small with current data (see right panel 
of Fig. I^for a comparison of these approximations in the 
parameter range of interest). 

Since the energy density of gravitational waves red- 
shifts as radiation (so long as we can neglect the contri¬ 
bution of the long wavelength modes) the relative energy 
density of gravitational waves compared to that of pho¬ 
tons, and thus the gravitational wave contribution to the 
value of iVeff, remains constant for all times of interest. 

In order to estimate the size of the effect of the ad¬ 
ditional radiation energy density, we will calculate the 
change in the redshift of matter-radiation equality. At 
the instant of matter-radiation equality, we have by def¬ 
inition 


Pm{f^eq) — Pr{^e<i) • 


(13) 


This can be written as an equation for Oeq as follows 


SttG 


eq 


4/3 


1+3(11) 



(14) 


where ag is the radiation constant, T^, 0 is the CMB tem¬ 
perature today, Aefi is given by Eq. (11), and we have 
fixed oq = 1. 

In order to place constraints on r and rit we have to 
make an assumption about the ratio fcuv • In our anal¬ 
ysis of the data we chose a pivot scale A:* = 0.01 Mpc“^. 
For LIGO we do not have to make an assumption on the 
UV cutoff, however we do assume the scaling is power law 
all the way up to fcLiGO- Th® minimum for this is ratio 
would then be fcuv/fc* — 10^°, while if we assume that 
the power law spectrum extends over ^ 60 e-folds^ the 
ratio is given by fcuv/fc* — 10^^. We estimate the con¬ 
tribution of primordial gravitational waves to the total 
energy density at the time of equality between matter 


^ This value is chosen not due to some particular mechanism for 
generating the gravitational waves in the early universe, but in¬ 
stead to roughly correspond to the maximum amount of hot big 
bang expansion. 


and the relativistic components of the universe. Shifts 
in the redshift of matter-radiation equality by roughly a 
percent will be constrained by the CMB power spectrum, 
since this the accuracy with which the redshift of equal¬ 
ity is measured. In Fig. we plot the relative change 
of the Ceq as a function of rit for these two numbers. It 
can be seen from this estimate that even for the lower 
value of fcuvi the observational constraint from the mod¬ 
ified expansion history will be at a level similar the best 
constraint from the CMB-I-LIGO data. 

For consistency we include this effect in our CMB anal¬ 
ysis. We add the gravitational wave contribution to the 
total the massless degrees of freedom in TVeff as computed 
above. The way the sampler usually scans the likelihood 
is to sample the densities and the optical depth inde¬ 
pendently (slow parameters) of the primordial param¬ 
eters (fast parameters); this scheme obviously does not 
work in case one of the densities depends explicitly on the 
primordial parameters. For that reason we changed the 
sampling. Unfortunately this makes the scanning of the 
likelihood slower, as one adds three additional parame¬ 
ters to the slow part of the code (namely the computation 
of the radiation transfer functions). In order to speed up 
the convergence, we therefore fix the foregrounds to their 
best fit values. 

The approximations we have used are valid for a pos¬ 
itive tensor tilt, and only for small values or negative 
values of the tilt would we have to worry about IR cor¬ 
rections. However, from the analytical analysis, the back¬ 
ground will not be affected in that limit and hence the 
total energy density from gravitational waves will be neg¬ 
ligible. The code does not compute a correction to Neff 
when rit <0. For very large values of rif, Eq. El breaks 
down predicting negative values of Neff, and so in that 
case, we replace A^eff with a very large number. The 
results are shown in Fig. [7] for two assumptions about 
the UV cutoff; one in which the cutoff is set at the 
LIGO scale and one in which we assume the maximal 
cutoff kfjy/kt — 10^^. The constraints are summarized 
in Tab. |lj Compared to our analytical estimate, we find 
the constraints to be slightly weaker. Although Zeq is con¬ 
strained quite strongly, Hq and the matter density tend 
to vary as rit increases in order to keep Zeq fixed, which 
therefore gives a weaker constraint than our naive pre¬ 
diction. This is shown in Fig. The strongest bound, 
as expected, is derived when we include BICEP2 data 
(which drives r to large values) and fcuv/^* — 10^^. We 
find nt = 0 . 331 ^: 25 . This stronger bound is driven by a 
non-zero r. Without the addition of BICEP2 data con¬ 
straints are weaker and are found to be remarkably close 
to LIGO bounds. 

The change in the expansion history due to the gravita¬ 
tional wave background energy density also affects BBN. 
The presence of extra massless degrees of freedom during 
BBN increases the primordial Helium fraction by increas¬ 
ing the expansion rate and thus decreasing the time dur¬ 
ing which free neutrons can decay before becoming bound 
into light nuclei |3T]. In m and recently in [55] it was 
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claimed that BBN puts an upper limit (with unknown 
confidence) on rit < 0.15. However, this constraint is de¬ 
rived under the assumption that fcuv/fc* — 10®^ > which 
would suggest the power law spectrum of gravitational 
waves extends over more than 120 e-folds and that many 
modes would have certainly been super-Planckian during 
the hot big bang expansion of the universe. Without a 
consistent theory of quantum gravity, it is unclear how 
to treat the evolution of such high energy gravitons. If 
one repeats the analysis with more sensible numbers, we 
find rit < 0.4 for r — 0.2, in line with and [21]. This 
bound again is very close to the bound from LIGO and 
from the CMB, and for completeness it should be taken 
into account. However, we expect based on a very sim¬ 
ilar analysis in [38] that the modified expansion history 
leads to effects for which the CMB provides a stronger 
constraint, given the same UV cutoff. That being said, 
the BBN constraint depends on the value of both Yp and 
fib and it would be possible to include this fully self con¬ 
sistently inside CAMB. In this analysis we have fixed the 
value oiYp = 0.24, and we leave a full treatment includ¬ 
ing these effects to future work. 

When considering the constraints derived from the 
contribution of gravitational waves to the massless de¬ 
grees of freedom, there are two important things to re¬ 
member. Firstly, the total energy density in gravitational 
waves depends on the entire primordial spectrum, which 
we have assumed here to be a power law all the way up 
to some UV cutoff. A more gradual transition to the re¬ 
heating phase, or production of additional gravitational 
waves during reheating could alter the constraints de¬ 
rived here. In addition, since the constraints on modified 
expansion history are from ‘indirect’ measures, the effect 
of gravitational waves would be degenerate with modi¬ 
fications to the neutrino sector or the addition of dark 
radiation which would have the same phenomenological 
effect on the measured CMB |5B]. As such, this con¬ 
straint is less robust than actually putting upper bounds 
on direct gravitational wave detection. In our analysis 
we have assumed that beyond photons and three families 
of standard model neutrinos, all of the radiation energy 
density is made up of gravitational waves, with no addi¬ 
tional massless degrees of freedom making a significant 
contribution to the energy density of the universe. 


IV. CONCLUSION 

In this paper we performed the first joint analysis of 
CMB data and late time measures of gravitational waves. 
This study is interesting for the following reasons. First, 
a combined analysis allowed us to put a bound on both 
the tilt and the running of the primordial tensor power 
spectrum. Secondly, even without BICEP2 data or the 
BICEP2/KECK and Planck cross data, we can put a 
bound on the tensor-to-scalar ratio r that is close to the 
constraint without allowing for the tensor tilt to vary. 
This is especially important given the recent conclusion 


that dust could be responsible for all of the detected B- 
mode power. While additional sources of gravitational 
waves could explain the detection of r (even without con¬ 
sidering dust), this gives interesting bounds on the energy 
scale of inflation [32] , of the detailed perturbation physics 
in this early epoch [9] , and of other early universe physics 
m- Thirdly, it has shown that all possible constraints 
on the tilt (and the running) are of the same order, and a 
full analysis should consider all these effects at the same 
time. 

Current data prefers the tilt to be positive and the 
running to be negative. The posterior parameters ranges 
in the combined analysis are reduced by factor 5 for the 
tilt and infinitely for the running (there is no constraint 
from the CMB). As expected, the data is not yet capa¬ 
ble of putting the inflation consistency condition to the 
test, but this work suffices as a confirmation that such an 
analysis is possible and that the combined data is already 
much more constraining than the CMB alone. A true test 
of the inflation consistency condition requires improved 
CMB data and measurement of gravitational waves at 
frequencies that are not dominated by foregrounds. 

We also showed that for large positive values of the 
tensor tilt, one expects the background evolution to be 
modified, altering the peak structure of the CMB. There 
are some caveats in this analysis, since it requires an 
explicit assumption about the UV cutoff of the tensor 
spectrum and one could in principle introduce other in¬ 
gredients with degenerate cosmological effects, i.e. this 
measure is not as ‘clean’ as direct constraints obtained 
with gravitational wave detectors. Eor that reason, one 
should be careful interpreting these results. Although an¬ 
alytical results indicated that the constraints would be 
better, performing the actual analysis showed that pa¬ 
rameter degeneracies lead to bounds on the tilt that are 
slightly worse (fcuv/fc* = 10^°) or only marginally better 
{k\jy/k^ = 10^^) than constraints from direct gravita¬ 
tional wave detectors. 

Constraints on nt from direct CMB measurements are 
obviously dependent on a non-zero detection of r, and 
are unlikely to greatly improve with future measurements 
given the large foregrounds and the limited lever arm for 
B-modes. On the other hand, direct measures place use¬ 
ful upper bounds on nt which are not particularly sensi¬ 
tive to the value of r. In the near future we expect LIGO 
to improve its bound on gravitational waves and as such 
on nt- At the same time, the background expansion will 
also be measured more accurately, with an expected error 
on the number of massless degrees of freedom of order a 
percent (see for example m)- Hence, we expect both 
measures to be equally important. In the long term how¬ 
ever, testing the consistency condition would require a 
more ambitious direct detection experiment with a much 
higher level of precision. For that reason alone, we antic¬ 
ipate the analysis presented in this paper to be of interest 
in the future, especially if the value of r is constrained at 
the percent level. 
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FIG. 7: The constraints on the gravitational wave spectrum when adding in the change in A'eff from gravitational waves. 
Allowing for the A'efi contribution cuts off the distribution on rit to more positive values. Adding in BICEP2 data again moves 
the constraints to more positive values of r, as expected. 'When considering the BKxP data in addition to Planck, the upper 
limits on nt are hardly affected as can been see in the hgure on the right. We also show the effect of changing the cut-off scale 
from kt/kyv — 10^° to fc*/fcuv = 10^"^. The overall effect is to limit the upper bound on nt- 
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Appendix A: On the gravitational wave energy 
density 

In this appendix we consider the definition of the en¬ 
ergy density of gravitational waves, following the treat¬ 
ment of [2211431 [521 EH. We are ultimately interested in 
how gravitational waves curve the background spacetime 
and impact the expansion rate. We begin by splitting 
the metric into a smooth background and a perturbation 

where we take the background to be given by the flat 
Friedmann-Robertson-Walker metric and the perturba¬ 
tion is defined in terms of the gauge invariant tensor per¬ 
turbation as 


5gij = , (A2) 

with hij^j = 0 and ha = 0. Understanding the impact 
of the gravitational waves on the evolution of the back¬ 
ground requires that we consider only fluctuations for 


which k/a{t) ^ H{i), modes which are deep inside the 
horizon. This might seem to pose a problem for the case 
we are considering, where there exist gravitational waves 
with wavelengths comparable to and larger than the hori¬ 
zon size at the epochs of interest. We will show below, 
however, that such modes give a negligible contribution 
to the energy density and can be safely neglected for the 
cases we consider. 

Now let us expand the Einstein equations around the 
background metric. There exist two small expansion pa¬ 
rameters; the first is the amplitude of the gravitational 
waves 5g = 0{\5g^^\), and the second is k/aH. The Ricci 
tensor can be expanded up to second order in Sg^iu as 

Rnf = + R^il} + R^^v + ■ • • ■ (A3) 


Now the key point is that R^lv is constructed purely from 
g^^i, and therefore only contains only low frequency modes 
(by which we mean variations on scales comparable to 
H~^), while is linear in Sg^iv and therefore contains 
only high frequency modes. However, is quadratic in 
Sg^ti, and therefore contains both high and low frequency 
modes. We can therefore split the Einstein equations into 
low and high frequency parts 


Rn-y — 




Low 


-b SttG ( - -g^ivT 


Low 


(A4) 


and 




High 


High 

+ 8^G ( ) . (A5) 


Eq. (A51 is a wave equation which governs the propaga¬ 
tion of 6g^i, on the background and is not of immediate 
interest for us. Eq. (A4| describes how the presence of 
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FIG. 8: A scatter plot showing the relation between the tensor tilt nt, the Hubble rate Hq and the redshift of equality Zeq, 
including the contribution from gravitational waves to A^efi- The top panel shows the case where the data are combined with 
BICEP2, while the middle and lower panels replace the BICEP2 data with the BICEP2/KECK-PlancA: cross spectrum for two 
different values of the cutoff. In all cases, «eq is fairly well constrained, hence in order to accommodate large nt (and hence a 
large A^eff, Hq is increased by adding more matter to the universe. Because of the steep increase of as a function of nt 
there is a very strong turn in degeneracy between nt and Hq. There is a shift in the upper bound of nt allowed by BICEP2 
compared to the cross spectrum, and this limit also shifts with different values of the cut off, as can be seen by comparing the 
middle to lower panels of the figure. 


gravitational waves curves spacetime and allows us to de¬ 
fine the energy momentum tensor for gravitational waves. 

Since we are considering only perturbations which vary 
on scales much smaller than the horizon size, we can 
project onto low frequency modes by performing an av¬ 
erage over a region which contains several wavelengths 
of the gravitational waves (or equivalently over a time 
in which the gravitational waves undergo many oscilla¬ 
tions). We therefore find 


This quantity can be computed explicitly in the trans¬ 
verse traceless gauge as 

+ 0{5g ^), (A8) 

where a vertical bar indicates a covariant derivative with 
respect to the background metric. The energy density of 
gravitational waves is then given by 


= -(4-) + 87rG(r^, - , (A6) 

and we can see that gravitational waves have an effective 
energy momentum tensor which is given by 

+ 0{5g^). (A7) 


= - ‘^H)5gt,{d^ - 2H)5gM) 

+ 0{5g^) 

= ■ (A9) 
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With the power spectrum A^(ry, k) defined as 

{hij{r],x)h’-^7],K)) = J dlogfcA^(? 7 ,fc), (AlO) 

and the transfer function T{r], k) defined as 

^h(k,ri) = Pt(k)[T(jn,k)f , (All) 

the energy density of gravitational waves can be written 
as 


PGW = I d\ogkPt{k) [r{v,k)f , (A12) 


where a prime refers to a derivative with respect to con¬ 
formal time rj. 

Constraints on direct detection of gravitational waves 
tend to be quoted in terms of the normalized energy den¬ 
sity per logarithmic scale 


^Gwik) 


1 dpGW 

Pcrit,o dlogfc 


Ptjk) 


[r{rio,k f , (A13) 


where we have used the fact that the critical energy den¬ 
sity today is given by Pcrit,o = SHq/SttG. 


Now let us briefly comment on the limits of the integral 
over wavenumber appearing in Eq. (A12). In order to be 
consistent, we should fix the lower limit to be at a scale 
for which fciR,/a(t) H{t) since we restricted ourselves 
in this derivation to modes which were deep inside the 
horizon. On the other hand, the derivative of the transfer 
function scales as in the limit that k —>■ 0 [23] . This be¬ 
havior ensures that the integral will converge in the limit 
km 0 for nt > —4. Furthermore the long wavelength 
modes (those for which k/a < ^H) will make a negligi¬ 
ble contribution to the total energy density for n* > —2. 
Since we are mostly interested in cases where rit > 0, 
we can safely take km — 0 without introducing a large 
error. The upper limit on the integral generally depends 
on the mechanism by which the gravitational waves were 
produced in the early universe. We consider a few cases 
in the main text and discuss how our constraints depend 
on this choice. 
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